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Determination of waveguide parameters 

Field of the Invention 

The present invention relates to the determination of parameters of a waveguide from 
5 wavefield data acquired from waves propagating in the waveguide. The invention may 
be used to obtain information about, for example, the thickness of the waveguide and/or 
the velocity of propagation of waves within the waveguide. 

Background of the Invention 

1 0 Guided waves result from rnultiple reflections in layered media. A guided wave occurs 
when a wave propagating in a layer is incident on a boundary of the layer at an angle to 
the surface normal greater than the "critical angle". As is well known, when the 
incident angle to the surface normal is greater than the critical angle, the sine of the 
angle of transmission to the surface normal, as determined by SnelTs Law, is greater 

15 than 1 . A wave transmitted out of the layer may exist only as an evanescent wave, and 
the phenomenon of "total internal reflection" occurs. The waves are therefore trapped 
in the layer and propagate as guided waves within the layer that do not decay 
significantly with travel distance. A layer that transmit guided waves in this way is 
called a "waveguide". The waveguide effect is shown schematically in Figure 1. A 

20 wave 1 propagating in a layer 2 is incident on a boundary 3 of the layer. When the 
incident angle to the surface normal is greater than the critical angle 9c, total internal 
reflection occurs leading to a guided wave. This is shown by the full line in Figure 1. 
(If the angle of incidence, to the normal to the surface 3, is less than the critical angle 
then the wave is partially transmitted and partially reflected, as shown by the dotted 

25 lines.) 

When the guided waves are recorded along the waveguide surface the apparent 
slowness or phase velocity of the guided waves exhibit dispersion - that is, they are 
each a function of the frequency of the waves. This dispersive character of guided 
30 waves can be used to obtain information about a waveguide or, more generally, about 
material properties of, and wave propagation velocities in, layered media. This is 
commonly achieved by solving an inverse problem, and matching a dispersion curve 
obtained for a numerical model of a waveguide to an observed dispersion curve. One 
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such method is described by M. Roth and K. Holliger, in "Joint inversion of Rayleigh 
and guided waves in high-resolution seismic data using a genetic algorithm", Society of 
Exploration Geophysicists, Expanded Abstracts, ppl 570-1 573 (1998). They suggested 
picking the dispersion curves of Rayleigh waves and guided waves in dispersion images 
5 and inverting them to obtain the velocity of P-waves, the velocity of S-waves, and the 
density of the waveguide using a genetic algorithm. This method requires an iterative 
inversion for P-velocity, S-velocity and waveguide density in the top layer and 
halfspace - 6 parameters in total - and so requires considerable computing power. 
Furthermore, this method requires the waveguide density, either as prior knowledge or 
10 as a parameter in the inversion. 

Dispersion images have also been used, by J. Xia et al, in "Estimation of near-surface 
shear-wave velocity by inversion of Rayleigh waves", Geophysics, Vol. 64, pp691-700 
(1999), to analyse dispersion curves of Rayleigh waves. The inversion of Rayleigh 
15 waves provides only S-velocities. The P-velocities are assumed to be known. This 
method also requires knowledge of the density of the waveguide layer. 

Summarv of the Invention 

20 In many technological applications of wave propagation, such as seismic surveying, one 
particular interest is the determination of wave propagation velocities. In many cases 
the earth's interior has a layer at or near the surface in which the velocity of propagation 
of seismic energy is different from the velocity of propagation of seismic energy in 
underlying layers. The presence of this layer causes a shift, known as the "static shift", 

25 in the arrival time of seismic energy recorded at seismic receivers disposed on or at the 
earth's surface, compared to the arrival time that would be recorded if the seismic 
velocity in the surface or near-surface layer were the same as the seismic velocity in the 
underlying layer. The static shift must be taken into account when analysing seismic 
data, and this requires knowledge of the velocity of propagation of seismic energy in the 

30 surface or near-surface layer and the thickness of the surface or near-surface layer. One 
approach to measuring these quantities is to treat the surface or near-surface layer as a 
waveguide and estimate the seismic velocity in the surface or near-surface layer and the 
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thickness of the surface or near-surface layer from measurements of the wavefield 
recorded at the earth's surface. 

The present invention provides a method of determining at least one parameter of a 
5 waveguide from wavefield data acquired from wave propagation in the waveguide, the 
method comprising the steps of: obtaining first and second dispersion curves in the 
frequency domain from the wavefield data; and determining at least one parameter of 
the waveguide from a frequency interval between the first dispersion curve and the 
second dispersion curve. 

10 

The invention uses the multiple character of guided waves to obtain waveguide 
parameters, such as a wave propagation velocity within the waveguide (hereinafter 
referred to as a "waveguide velocity" for convenience") and/or the thickness of a 
waveguide, from measurements of the wavefield recorded at the surface of the 

15 waveguide. If the waveguide thickness is of the order of the wavelength or thinner, the 
multiple reflections of the guided waves interfere with one another and so are not 
separable in the time domain. The invention overcomes this by making use of the fact 
that, in the fi*equency domain, the dispersion curves corresponding to different guided 
wave modes are separate. Parameters of the waveguide, such as waveguide velocity 

20 and the thickness of the waveguide may be determined from the dispersion curves in the 
fi*equency domain. 

The invention can be applied to any wavefield that propagates in a waveguide and is 
recorded at the surface of the waveguide. In the field of seismic surveying, for example, 

25 it may be used to estimate parameters of an overlying layer of the earth's interior for use 
in correction of the static shift, or in wavefield separation. In a land-based seismic 
survey, an elastic wavefield is generated by explosive devices or vibrators and the 
wavefield propagating through a near surface layer that acts as a waveguide is recorded 
by geophones. In marine seismic surveying, an acoustic wavefield is commonly 

30 generated in the water by airguns and recorded in the water by hydrophones. If the 
wavefield is recorded on the seabed, hydrophones and geophones are used to record the 
acoustic and elastic wavefields respectively. 
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When applied to seismic surveying, the invention is able to estimate the velocity of P- 
waves and/or the velocity of S-waves. The method of the invention is independent of 
the density of the medium, and so does not require this information to be known. 

5 A second aspect of the invention provides a method of processing wavefield data, the 
method comprising: acquiring wavefield data; determining at least one parameter of a 
waveguide according to a method of the first aspect of the invention; and taking the at 
least one parameter into account during subsequent processing of the wavefield data. 

1 0 The wavefield data may be seismic wavefield data. 

A third aspect of the invention provides an apparatus for determining at least one 
parameter of a waveguide from wavefield data acquired from wave propagation in the 
waveguide, the apparatus comprising: means for obtaining first and second dispersion 
1 5 curves in the frequency domain from the wavefield data; and means for determining at 
least one parameter of the waveguide from a frequency interval between the first 
dispersion curve and the second dispersion curve. 

The apparatus may comprise a programmable data processor. 

20 

A fourth aspect of the invention provides a storage medium containing a program for 
the data processor of an apparatus as defined above. 

A fifth aspect of the invention provides a storage medium containing a program for 
25 controlling a programmable data processor to carry out a method of the first aspect of 
the invention. 

Brief Description of the Drawings 

Preferred embodiments of the present invention will now be described by way of 
30 illustrative example with reference to the accompanying figures, in which: 

Figure 1 is a schematic illustration of propagation of a guided wave; 
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Figure 2 is a schematic illustration of propagation of seismic energy in a layer of the 
earth's interior; 

Figure 3(a) shows seismic data acquired in the seismic surveying arrangement of figure 
5 2; 

Figures 3(b) and 3(c) illustrate steps in obtaining a dispersion curve in the frequency 
domain from the seismic data of figure 3(a); 

10 Figure 4 illustrates a typical set of dispersion curves obtained fi*om the seismic data of 
figure 3(a); 

Figure 5 illustrates one method of determination of waveguide parameters from the 
dispersion curves of figure 4; 

15 

Figures 6(a) to 6(f) illustrate the present invention compared with processing in the x-V 
domain; 

Figure 7(a) illustrates a typical dispersion image in the fi*equency domain; 

20 

Figure 7(b) illustrates the auto-correlation of the image of figure 7(a); 

Figure 7(c) to 7(i) illustrate converted images of the auto-correlation image of figure 
7(b) for different assumed layer velocities; 

25 

Figure 8 is a schematic flow diagram of one embodiment of the present invention; 

Figure 9 is a schematic flow diagram of a second embodiment of the present invention; 

30 Figure 10 illustrates the second embodiment of the present invention; 

Figure 1 1 is a schematic flow diagram of a third embodiment of the present invention; 
and 

5 
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Figure 12 is a block schematic diagram of an apparatus according to the present 
invention. 

5 Detailed Description of the Invention 

Figure 2 illustrates typical paths of seismic energy at a survey location where the earth's 
interior is represented by a layer 5 that is overlaid by a thin layer 4. A seismic source 6 
is located in the overlying layer 4 and, when actuated, emits seismic energy. The 
emitted seismic energy is detected by an array of suitable detectors 7 disposed at the 
10 earth's surface. 

The overlying layer 4 has thickness acoustic velocity cj and density p/. The 
underlying layer 5 has an acoustic velocity C2 and density p2, which are assumed to be 
different from the acoustic velocity and density of the overlying layer 4 (c^ will be 
15 greater than cy). In theory, the layer 5 is assumed to have infinite depth. However, in 
practice, described methods also work when layer 5 has finite depth. 

Figure 2 shows three possible paths of seismic energy emitted by the source 6. Seismic 
energy that travels along the path shown as a dotted line is incident on the interface 8 
20 between the overlying layer 4 and the underlying layer 5 at an angle, to the normal to 
the interface, that is less than the critical angle. As a consequence, some of the seismic 
energy is reflected and some is transmitted. 

Seismic energy that travels along the path shown as a broken line is incident on the 
25 interface 8 at exactly the critical angle to the normal to the interface 8. This energy 

undergoes critical refraction, propagates along the interface 8, and is subsequently 
refracted upwards. 

Seismic energy that travels along the path shown as a full line is incident on the 
30 interface 8 at an angle to the normal to the interface that is greater than the critical 
angle, and undergoes total internal reflection. This energy is trapped in the layer 4, and 
forms a guided wave that propagates through the layer 4. The layer 4 thus acts as a 
waveguide. 
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Figure 3(a) shows the seismic energy recorded at the receivers 7 in figure 2. Figure 3(a) 
shows the amphtude of the vertical velocity component of the acquired seismic energy 
as a function of the horizontal distance X from the source and as a function of the time t 
5 since the actuation of the source 6. The greatest positive amplitude in figure 3(a) is 
shown as white, and the greatest negative amplitude is shown as black, as indicated in 
the key at the right of figure 3(a). It will be seen that events corresponding to reflected 
waves (such as the path shown in dotted lines in figure 2), to critically refracted waves 
(such as the path shown in broken lines in figure 2) and to guided waves (such as the 
10 path shown in full lines in figure 2) can be identified in the seismic data. Figure 3(a) is 
reproduced as an inset in figure 2, and the dotted arrow, the broken arrow and the full 
arrow in the inset in figure 2 indicate events corresponding to, respectively, reflected 
waves, critically refracted waves and guided waves. 

15 The seismic data of figure 3(a) may be recorded in any conventional manner. Typically, 
the data are obtained using digital receivers that sample the wavefield at the receiver 
position at regular intervals. The sampled data may be stored in each receiver for 
subsequent retrieval, or they may be transmitted immediately to a central computer for 
storage there. 

20 

The recorded data of figure 3(a) are the vertical component of the velocity field as a 
function of the x-coordinate and time, and may be denoted as s(x,t). In order to perform 
the method of the invention, the recorded data s(x,t) are transformed into a dispersion 
image in the fi*equency domain. Figure 3(c) shows the result of transforming the data of 

25 figure 3(a) into a dispersion image in the frequency-phase velocity domain via the t-V 
domain as shown by McMechan and Yedlin in "Analysis of Dispersive Waves by 
Wavefield Transformation" Geophysics Vol. 48 pp869-874 (1981). In figure 3(c), the 
regions of highest amplitude in the dispersion image are shown in black, while regions 
of lowest amplitude are shown as white, as indicated by the key at the right of figure 

30 3(c). It will be seen that the regions of highest amplitude in the dispersion image of 
figure 3(c) align along dispersion curves, and these are indicated in figure 3(c) as full 
lines. Each dispersion curve 9a, 9b, 9c shown in figure 3(c) corresponds to a different 
mode of wave propagation. 
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The curves shown as full lines in figure 3(c) are calculated according to equation (1) 
below. 

Figure 3(b) illustrates, as a comparison, the data in the time domain, in this case in the x- 
V (intercept time-velocity) domain. It will be seen that the different reflection events 
interfere, and are not separable fi*om one another. 

Figure 3(b) represents an intermediate step in the transformation to the frequency 
domain. The image in the frequency-velocity domain is, in this example, obtained by 
the method of McMechan and Yedlin (above), who first transform data to the t-p 
(intercept time-slowness) domain, and then apply a Fourier transform to get to the f-v or 
f-p domain. An alternative approach is to transform the data from the x-t domain to the 
f-x domain, and then to the f-p or f-v domain. 

A dispersion relation for acoustic guided waves has been derived by W M Ewing et al in 
"Elastic Waves in Layered Media", McGraw-Hill (1957). They derived the following 
expression for the dispersion curves of acoustic-guided waves: 




In Equation (1), V is the phase velocity, and (o (= 2;z/) is the angular frequency: The 
parameters h, cj, pj, C2 and p2 have the meanings defined with reference to Figure 2. 

Equation (1) correctly predicts that, at a given phase velocity, the frequency interval 
between each pair of adjacent or successive dispersion curves is constant since the 
argument of the tangent function must fulfil: 
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Equation (2) indicates that the frequency difference between successive dispersion 
curves is a function of the phase velocity and is given by: 



The number n refers to the mode number in dispersion curve theory. 

It will be noted that the densities pi, p2 do not appear in equation (3). The present 
invention takes advantage of this, and uses equation (3), rather than equation (1), for 
waveguide parameter estimation. 

Figure 4 shows a set of dispersion curves having the general formula of equation (1). It 
was shown in equation (3), that, at a given phase velocity V, there is a constant 
frequency difference between each pair of immediately successive dispersion curves. It 
will be seen that the dispersion curves of figure 3(c), obtained from the seismic data of 
Figure 3(a), have the same form as the theoretical dispersion curves shown in figure 4. 

The present invention provides a method for obtaining properties of a waveguide from 
dispersion curves that have the general form shown in figure 3(c) or figure 4. The 
essential concept of the invention is to determine the frequency difference 4/^ between 
adjacent dispersion curves (i.e., between dispersion curves of adjacent modes) for a 
particular phase velocity. It is then possible to determine parameters of the waveguide 
in which the wavefield data used to obtain the dispersion curves were obtained. 

In a simple realisation of the invention, the frequency difference zl/^ between adjacent 
dispersion curves may be simply picked from a dispersion image that is displayed, for 
example, as a hard copy or on a computer screen. 




(3) 
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In order to solve equation (3) to determine the thickness h of the waveguide expHcitly, it 
is necessary to know the frequency difference Af for a particular phase velocity, and to 
have an estimate of the waveguide velocity cy. Both these quantities may be determined 
5 from a dispersion curve of the type shown in figure 4. As is clear from figure 4, the 
waveguide velocity cj is given by the asymptotic velocity of the dispersion curve(s) at 
high frequency values. The asymptotic velocity of the dispersion curve is shown as a 
broken line in Figure 4. 

10 This embodiment of the invention is illustrated in figures 5 and 8. Figure 5 illustrates 
how the quantities Af (for a given phase velocity) and cj may be determined from 
displayed dispersion curves. In principle, the frequency difference Af between two 
successive dispersion curves is simply measured from displayed dispersion curves for a 
known value of the phase velocity. If only two dispersion curves are available, then a 

15 single value for the frequency difference at a particular phase velocity can be obtained. 
However, if three or more dispersion curves are available, it is possible to obtain two or 
more values of the frequency difference for a known value of the phase velocity. In 
figure 5, for example, three dispersion curves are visible, which have frequencies at a 
particular phase velocity V of fufi and fs respectively. It is therefore possible to obtain 

20 two values for the frequency difference at a given phase velocity F, as follows: 

4/i=/2-/7; 

25 These two values of the frequency difference at a particular phase velocity may be 
averaged using 4/av (4/i *^ 4^)/2. 

As noted above, the asymptotic velocity to the dispersion curves is equal to the 
waveguide velocity c/ and this may, in principle, again be directly determined from a 
30 displayed dispersion curve as shown in figure 5. Thus, the value of the waveguide 
velocity may be directly obtained from a dispersion curve in the frequency domain (or 
may be determined as the average of the asymptotic velocities of two or more dispersion 

10 
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curves). Furthermore, once the waveguide velocity c/ and the frequency difference 
Af{V) between successive dispersion curves at a known phase velocity V has been 
determined, it is then possible to determine the waveguide thickness from equation (3), 
All quantities in equation (3), except for h, may be obtained from the dispersion image. 
5 In this connection, it will be noted that the density of the waveguide does not appear in 
equation (3). 

Figure 8 is a flow diagram showing the principal steps of the method illustrated in 
figure 5. Initially, at step 81, a wavefield is recorded in the time-space domain. This 
10 may be done, for example using an array of receivers as shown in figure 2 disposed at 
the surface of a waveguide. This step will typically record the amplitude of a 
component of the wavefield as a function of space and time. Figure 3(a) illustrates one 
form in which these results may be obtained. 

15 Next, at step 82, the wavefield data obtained in step 81 are transformed from the space- 
time domain into the frequency domain. Step 82 may consist of transforming the data 
into the frequency and phase velocity domain as in figure 3(c), or into the fi*equency and 
slowness domain. The transform may be carried out using any suitable technique. The 
result of step 82 is a dispersion image that will be generally similar to the dispersion 

20 image of figure 3(c). A set of dispersion curves may then be defined, for example as 
curves defined by the regions in which the amplitude of the dispersion image is a 
maximum. 

Once a set of dispersion curves has been obtained, the asymptote to the dispersion 
25 curves is measured. If step 82 provides a dispersion image in the fi-equency and phase- 
velocity domain, step 83 will provide an asymptotic value of the phase velocity which, 
as explained above, is equal to the phase velocity of the waveguide. If step 82 provides 
a frequency-slowness dispersion image, the asymptote corresponds to the slowness of 
the waveguide. There is a corresponding version of equation (3) for this case, in which 
30 V is replaced by px and c/ is replaced by i/p, where p is the medium slowness. 
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Next, at step 84, the frequency distance Af between two successive dispersion curves is 
measured. This is measured for a known value of the phase velocity if step 82 provides 
a dispersion curve in the frequency and phase velocity domain, or is measured at a 
known value of the slowness if step 82 provides a dispersion curve in the frequency and 
5 slowness domain. If more than two dispersion curves are available, more than one 
value of Af may be determined at step 84 and in this case an average frequency 
difference, Afav, can be determined. 

It should be noted that steps 83 and 84 do not need to be performed in the order shown 
10 in figure 8. These steps could be performed in the reverse order, and in principle could 
be performed simultaneously. 

Next, at step 85 the thickness h of the waveguide layer is determined using equation (3) 
using the value for the waveguide velocity cj determined at step 83, using the frequency 
1 5 difference Af (or average frequency difference determined at step 84, and using the 
value of the phase velocity V for which the frequency difference (or average frequency 
difference) was determined at step 84. 

Finally, at step 86 the value of the waveguide thickness h determined at step 85 and the 
20 value of the waveguide velocity cj determined at step 83 may be output. 

It should be noted that the wave velocity C2 of the layer below the waveguide layer (for 
example the layer 5 in figure 2) may also be determined from a dispersion curve of the 
general form of figure 3(c). As indicated in figure 4, this velocity is the maximum 
25 velocity at which the dispersion curve(s) exist. 

Of course, if prior information about the waveguide velocity cj is available, for example 
from knowledge of the composition of the waveguide, then this information may be 
used as well as, or instead of, information about the waveguide velocity cj derived from 
30 the dispersion curve. 
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The frequency difference is preferably determined at phase velocity where the 
resolution of the dispersion image is good. In principle, however, the frequency 
difference may be determined at any phase velocity between ci and C2. 



The simple method described in figures 5 and 8 can in principle provide good results if 
the dispersion image has a high resolution so that the frequency difference between 
adjacent dispersion curves and the asymptotic velocity can be determined accurately. 
However, this may not always be the case and, in particular, the asymptote to the 
dispersion curves may not be well defined so that the waveguide velocity c/ cannot be 
easily obtained by inspection of the dispersion curves. In such circumstances, it may be 
preferable to employ an alternative method that is shown in figure 9. 

In the method of figure 9, a wavefield is recorded in the space-time domain at the 
surface of a waveguide, and at step 92 this is transformed to give a dispersion image in 
the frequency domain. Steps 91 and 92 correspond to steps 81 and 82 of figure 8, and 
will not be described in further detail. 

Next, the fi-equency difference between two successive dispersion curve is measured for 
two or more different values of the phase velocity or slowness (depending on the 
particular dispersion image produced at step 92). For example, as shown in Figure 10, 
the frequency difference between two successive dispersion curves is measured at two 
different values Ky, F? of the phase velocity. If more than two dispersion curves are 
available, the frequency difference estimated for one or more of the phase velocities can 
be obtained using the averaging process described above. 

The thickness of the waveguide is then determined from the sets of values (Af (Ki), F]), 
(Af (K2), V2) etc., according to the following equation: 




(4) 
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Equation (4) essentially states that the same value of the thickness of the waveguide 
should be obtained from every set of (Af(Vi), K/), since the thickness of the waveguide is 
independent of the phase velocity. Equation (4) expresses this as a consistency criteria. 
If h and cj are both unknown, it is possible to determine both h and C/ provided that 
5 values for Af(V) are available for two or more different values of the phase velocity. 

It is very easy to obtain self-consistent values for the waveguide thickness and the 
waveguide velocity using equation (4), provided that the frequency difference is known 
for two or more different values of the phase velocity. In this connection, it will again 
10 be noted that the density of the waveguide does not appear in equation (4). It is 
therefore straightforward to obtain accurate results for the waveguide thickness and the 
waveguide velocity cy, even if the asymptote to the dispersion curves is not well-defined 
and the waveguide velocity cannot be determined by inspection of the dispersion 
curves. 

15 

In the method of figure 9, therefore, once values for the frequency separation between 
successive dispersion curves have been obtained for at least two different phase 
velocities, equation (4) is solved at step 94 to obtain self-consistent values for the 
waveguide thickness and the waveguide velocity. The self-consistent values of the 
20 waveguide thickness h and the waveguide velocity cj are those values that provide a 
value for the thickness h of the waveguide that is most nearly independent of the phase 
velocity. 

Figure 1 1 illustrates a further method of determining one or more waveguide parameters 
25 from wavefield data recorded in the space/time domain. Initially, wavefield data is 
recorded in the space/time domain at step 1 1 1 and is transformed to produce a 
dispersion image in the frequency domain at step 1 12. These steps correspond to steps 
81 and 82 respectively and will not be further described here. 

30 The method of figure 1 1 uses an auto-correlation technique to determine the waveguide 
thickness and the waveguide velocity. It is known to detect a periodic signal by 
calculating the autocorrelation of a series of samples and picking the maxima in the 



14 



57.0541 USNP 



auto-correlation function, but the invention extends this technique by applying it to the 
dispersion image obtained in the frequency domain. 

In this method the recorded data s(xj) in the time domain are transformed into a 
5 dispersion image in the frequency domain. For the purpose of description it will be 
assumed that the data are transformed into a dispersion image S(f,V) in the frequency- 
phase velocity domain, but the invention is not limited to this particular transformation. 
This transformation may be carried out according to any known technique such as, for 
example, the technique described by G A McMechan and M J Yedlin (above) or by C B 
10 Park et al. in "Imaging Dispersion Curves of Surface Waves on Multi-channel Record", 
Society of Exploration Geophysicists, Expanded abstracts, pp 1377-1380 (1998). This 
transformation produces complex-valued Fourier coefficients S(f,V) of the dispersion 
image. These complex-valued Fourier coefficients are then auto-correlated according 
to: 

15 

f^k,v ^^^iy^ky (5) 

This auto-correlation technique is adapted from the technique described by J F 
Claerbout in "Fundamentals of Geophysical Data Processing", Blackwell Scientific 

20 Publications (1976). In this method, the auto-correlation is applied in the frequency 
domain with frequency index / and frequency-lag index k for each value of the phase 
velocity V, The maxima of the modulus of the complex-valued auto-correlation 
function are located at integer multiples of Af, This is illustrated in figures 6(b) and 
6(d). Figure 6(b) shows the dispersion image obtained in the frequency and phase- 

25 velocity domain, and corresponds essentially to the image shown in figure 5, The 
results of the auto-correlation of the dispersion image of figure 6(b) are shown in figure 
6(d), with the maximum amplitude of the auto-correlation function being represented in 
black and with the minimum amplitude being represented in white as shown by the key 
at the right hand side of figure 6(d). The periodic structure of the results of the auto- 

30 correlation is apparent from figure 6(d) and, in principle, the frequency difference 
between adjacent branches of the auto-correlation image can be determined, as a 
function of the phase velocity, direct from the auto-correlated image. 
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Figure 6(a) and 6(c) correspond to figure 6(b) and 6(d) respectively, but show the t-V 
transformed data and its auto-correlation. It will be seen that auto-correlation in the 
time-domain cannot resolve the multiple structure of the dispersion image, since the 
5 multiple reflections interfere. In contrast, the multiple reflection modes are clearly 
separate in the frequency domain, as is clear from Figures 6(b) and 6(d), and the 
frequency difference between successive modes can be clearly measured in the auto- 
correlated image in the frequency domain. 

10 Figure 6(f) shows the results of the auto-correlation at a value for the phase velocity of 
1500m/s. It will be seen that there is a clear peak in Af at around 75-80 Hz and a 
smaller peak in Af at just over 150 Hz. This indicates that the value of Af is 
approximately 75 Hz at a phase velocity V= 1500 m/s. (The small amplitude peak at Af 
= 150 Hz is the result of correlation between the first and third branches or curves of 

15 the dispersion image shown in figure 6(b). 

Figure 6(e) shows the results of the auto-correlation in the x-V domain at a phase 
velocity of 1 500 m/s. It will be seen that there is no clear peak in this, and this is 
because of the overlap between the multiple reflections in the time domain. 

20 

Although the frequency difference for a particular phase velocity can be determined 
direct from the results of the auto-correlation in the frequency domain, the method 
shown in figure 1 1 uses a consistency criterion that is analogous to that used in the 
method of figure 9. In the method of figure 1 1, the consistency criterion is applied to 
25 the auto-correlation image of figure 6(d) by mapping each row of the matrix R^v from 

k - Af Xo h using equation (4). The auto-correlation event is aligned when the correct 
waveguide velocity cy is chosen. This process is illustrated in figure 7(a) to 7(i). 

Figure 7(a) and 7(b) correspond to figures 6(b) and 6(d) respectively, and show the 
30 dispersion image in the frequency domain (figure 7(a)) and the auto-correlation thereof 
(figure 7(b)). Figures 7(c) to 7(i) show the results of converting the auto-correlation 
image into the A-K domain (the waveguide thickness and phase velocity domain) using 
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equation (4) for seven different values of the waveguide velocity c/. Figure 7(c), for 
example, shows the results obtained by converting the auto-correlation image into the 
/z-F domain for an assumed waveguide velocity of Ci = 900 m/s, figure 7(d) shows the 
results obtained for an assumed waveguide velocity of ci = 1,000 m/s, and so on. 

5 

In figures 7(c) to 7(i) the results have been inverted to negative values of the layer 
thickness, so that the order of the events is the same as in the auto-correlation image. 

By inspection of figure 7(c) to 7(i), it is seen that the auto-correlation event is correctly 
10 aligned for a waveguide velocity of 1,100 m/s in figure 7(e). The auto-correlation 
image is not correctly aligned in the results obtained for other assumed values of the 
waveguide velocity. 

Furthermore, in figure 7(e) the auto-correlation event is aligned at a waveguide 
15 thickness h of 10m. Thus, the results of figure 7(c) to 7(i) indicate that the waveguide 
has a waveguide velocity cy of 1 ,100 m/s and a thickness h of 10m. 

The step of transforming the auto-correlation image to the thickness-phase velocity 
domain, and picking the correct waveguide velocity and thickness is step 114 in the 
20 method of figure 1 1 . The obtained values for the waveguide thickness and the 
waveguide velocity are then output at step 115. (In the case of figures 7(a) to 7(i), for 
example, the values h = 1 Om and cj = 1 ,1 00m would be output at step 1 1 5.) 

Where the invention is applied to a guided wavefield propagating in a layer at or near 
25 the earth's surface, such as the layer 4 of Figure 2, the invention is able to determine the 
wave propagation velocity in, and the thickness of, the surface or near-surface layer. 
These values may be used in subsequent processing of the wavefield. Where the 
recorded wavefield data are seismic data, for example, the invention may be used to 
determine one or more parameters of a layer of the earth, such as a surface or near- 
30 surface layer, that acts as a waveguide. These parameters may then be used in 
subsequent processing of the seismic data. For example, the waveguide parameters can 
be used to determine the static shift produced by the surface or near surface layer, and 
so allow the acquired seismic data to be corrected for the static shift. The waveguide 
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parameters may alternatively be used in wavefield separation. The parameters of the 
surface or near-surface layer may also be used in the processing of other seismic data 
acquired at that survey location. 

5 The invention has been explained above with reference to an acoustic waveguide. The 
invention is not limited to this, however, and may also be applied to a waveguide in an 
elastic medium. The dispersion curves of guided P-waves and guided S-waves can be 
analysed using the methods of the invention to obtain the P-wave and S-wave velocities 
of an elastic waveguide, and also to obtain the thickness of an elastic waveguide. The 
10 P-velocity of a waveguide can be determined using the above methods. The S-velocity 
of the waveguide can be determined from the dispersion image of the higher modes of 
the ground roll. The S-velocity may be determined using any of the methods described 
above, but applied to a different part of the data. To do this, the dispersion images must 
cover the range of P- and S- velocities, either in a single image or in separate images. 

15 

Any technique that creates a dispersion image from recorded wavefield data can be used 
to obtain the dispersion image in the frequency domain. One approach is to transform 
the recorded wavefield from the time-space domain into the frequency-slowness or 
frequency-phase velocity domain, for example as described by McMechen and Yedlin, 
20 above. A parametric spectral analysis tool such as the "FK-MUSIC" algorithm 
proposed by K Iranpour et al. in "Local Velocity Analysis by Parametric Wave Number 
Estimation in Seismic (FK-MUSIC)", EAGE Expanded Abstracts (2002) may be used. 
This technique can obtain a high-resolution dispersion image from fewer traces than 
with the method of McMechen and Yedlin. 

25 

The discussion given above has assumed a one-dimensional model of the waveguide, in 
which the waveguide is laterally invariant in velocity and thickness. Laterally-varying 
waveguide properties can be estimated using the methods of the present invention if the 
obtained dispersion image corresponds to a locally one-dimensional medium. Local 
30 dispersion images may be obtained with the FK-MUSIC algorithm of Iranpour et al. 

Figure 12 is a schematic block diagram of a programmable apparatus 10 according to 
the present invention. The apparatus comprises a programmable data processor 1 1 with 
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a program memory 12, for instance in the form of a read only memory ROM, storing a 
program for controlling the data processor 1 1 to perform any of the processing methods 
described above. The apparatus further comprises non-volatile read/write memory 1 3 
for storing, for example, any data which must be retained in the absence of power 
supply. A "working" or "scratchpad" memory for the data processor is provided by a 
random access memory (RAM) 14. An input interface 15 is provided, for instance for 
receiving commands and data. An output interface 16 is provided, for instance for 
displaying information relating to the progress and result of the method. Seismic data 
for processing may be supplied via the input interface 15, or may alternatively be 
retrieved from a machine-readable data store 17. 

The program for operating the system and for performing the method described 
hereinbefore is stored in the program memory 12, which may be embodied as a semi- 
conductor memory, for instance of the well-known ROM type. However, the program 
may be stored in any other suitable storage medium, such as magnetic data carrier 1 2a 
(such as a "floppy disc") or CD-ROM 12b. 

While the invention has been described in conjunction with the exemplary embodiments 
described above, many equivalent modifications and variations will be apparent to those 
skilled in the art when given this disclosure. Accordingly, the exemplary embodiments 
of the invention set forth above are considered to be illustrative and not limiting. 
Various changes to the described embodiments may be made without departing from the 
spirit and scope of the invention. 
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